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Stochastic  mixing  model  with  power  law  decay  of 

variance 

By  S.  Fedotov,  M,  Dime  and  H.  Pitsch 


1.  Motivation  and  objectives 

The  mixing  of  a conserved  scalar  c = c(t,x),  advected  by  a turbulent  flow,  remains 
a problem  of  both  fundamental  and  practical  interest.  One  of  the  basic  characteristics 
of  the  mixing  process  is  the  rate  at  which  the  scalar  variance  a^(t)  = ((c  — fi)^)  decays 
with  time.  Here  r,  is  the  mean  value,  and  the  angular  brackets  ( • ) denote  an  averaging 
procedure.  One  of  the  simplest  and  widely  used  mixing  models  is  the  interaction  by 
exchange  with  the  mean  (lEM)  (Villermaux  & Devillon  1972)  or  linear  mean  square 
estimate  model  (LMSE)  (Dopazo  & O’Brien  1974,  Sabel’nikov  & Gorokhovski  2001).  In 
this  model,  the  scalar  relaxes  toward  its  mean  r according  to  the  equation 


dc 

dt 


-i(c-/i). 


(1.1) 


In  the  coalescence-dispersion  model  (CD)  (Curl  1963),  mixing  of  two  particle  volumes 
is  described  by  (1)  a coalescence  and  (2)  dispersion  process  during  one  time  step  At. 
After  that  time,  both  particles  have  the  same  scalar  value,  which  is  equal  to  the  mean 
of  the  two  values  before  the  mixing.  The  pdf  relaxes  to  a bell-shaped  distribution  but 
deviates  from  the  Gaussian  shape  primarily  in  the  tails  of  the  distribution  (Peters  2000). 
Both  models  introduce  a time  scale  r which  is  commonly  approximated  by  a turbulent 
time  scale.  Other  mixing  models  are  based  on  mapping  closure,  Fokker-Planck  model, 
euclidean  minimum  spanning  trees  (EMST),  or  Langevin  model.  Reviews  of  these  can 
be  found  in  Peters  (2000)  and  Pope  (2000). 

There  are  two  different  laws  governing  the  decay  rate  of  a passive  scalar: 

(i)  the  exponential  law 

<r^(t)  oc  exp{-t/r}  (1.2) 

with  the  characteristic  time  r and 

(ii)  the  power  law 

<Tc"(t)oct-“  (1.3) 

without  any  characteristic  time  scale  (Lesieur  1997). 

Most  theoretical  models  introduce  a characteristic  time  scale  and  assume  implicitly  or 
explicitly  the  exponential  decay  rate  (1.2)  which  is  only  appropriate  for  stationary  tur- 
bulence. 

The  main  purpose  of  this  paper  is  to  study  the  mixing  process  following  the  power  law 
(1.3)  that  is  typical  for  decaying  turbulence.  Experimental  results  show  that  the  decay  of 
the  variance  strongly  depends  on  the  initial  ratio  of  the  velocity  and  scalar  length-scales 
and  that  there  is  no  universal  decay  exponent  (Durbin  1982). 

It  is  well  known  that  the  decay  exponent  a depends  on  the  low  wavenumber  part  of 
the  scalar  and  velocity  spectrum,  Ec{k,t)  and  E{k,t)  (Lesieur  1997).  One  can  expand 
both  spectra  into  a Taylor  series:  Ec{k,t)  = 27rfc^(Co  + <^2^^  -I-  0{k'^))  and  E{k,t)  = 
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2TTk^{Bo  + ^2^^  + (Chasnov  1994).  For  high  Reynolds  and  Peclet  numbers,  the 

scalar  variance  can  be  considered  as  a function  of  Co,Bq  and  t alone  and  dimensional 
arguments  lead  to 


alit)  oc 


(1.4) 


with  a = 1.2.  For  Co,  Bo  equal  to  zero,  the  scaling  argument  results  in  different  approx- 
imate decay  laws  with  a equals  to  6/7, 10/7  or  2 (Lesieur  1997). 

Here  we  present  a simple  stochastic  mixing  model  based  on  the  law  of  large  numbers 
(LLN)  (Feller  1966).  The  reason  why  the  LLN  is  involved  in  our  formulation  of  the 
mixing  problem  is  that  the  random  conserved  scalar  c = c(t,x{t))  appears  to  behave 
as  a sample  mean.  It  converges  to  the  mean  value  fj,,  while  the  variance  cr^(t)  decays 
approximately  as  Since  the  variance  of  the  scalar  decays  faster  than  a sample  mean 
(typically  a is  greater  than  unity),  we  will  introduce  some  non-linear  modifications  into 
the  corresponding  pdf-equation  (see  Eq.  (2.28)  below).  The  main  idea  is  to  develop  a 
robust  model  which  is  independent  from  restrictive  assumptions  about  the  shape  of  the 
pdf.  Here  we  exploit  the  similarity  of  the  behavior  of  a scalar  c to  that  of  the  sample 
mean 


Cn  — 


1 ^ 


k=l 


(1.5) 


where  Ci  i • • • i Cn  is  a sequence  of  mutually  independent  random  variables,  each  having  a 
mean  /x  and  standard  deviation  <t|.  The  LLN  tells  us  that  the  random  sum  Cn  tends  to 
the  mean  value  n with  probability  one,  while  the  variance  = {{cn  - nY)  decays  as 
cr|/n.  Then,  any  arbitrary  initial  pdf  p(no,c)  tends  to  a (J-distribution  as  n — ♦ oo. 

In  the  present  paper  the  discrete  increment  n can  be  understood  as  the  time  variable  t. 

The  main  result  of  this  paper  is  the  derivation  of  the  time-discrete  non-linear  integral 
equation  for  the  pdf  of  c 


1 

p{t  + l,c)  = j 
0 


i^_l 

f ^ + 1 (.  (t  -b  l)(c  - /i)  - A(1  + £)(??-  p) 

J 1 -(- 1 — A(1 -l- £■)  \ 1 -f  t — A(1 -|- s) 

-1 


p{t,Tf)ip{t,e)  de  dt] , 


(1.6) 


for  t = 1,2,3,...  . Here  is  the  pdf  for  the  exchange  rate  and  A is  the  mixing 

intensity  (see  below).  The  main  property  of  this  equation  is  that  it  describes  the  relaxation 
from  an  arbitrary  initial  distribution  to  a (5-function 


p{t,c) 6{c  — p.)  as  t — » oo  (1.7) 

and  the  decay  of  the  variance,  (7^{t)  = — p)^),  is  of  the  form  t"“.  The  case  a = 1 

corresponds  to  the  law  of  large  numbers. 

The  remainder  of  this  paper  is  organized  as  follows.  In  Section  2 we  derive  the  integral 
equation  from  a stochastic  difference  equation  describing  the  evolution  of  the  pdf  of  a 
passive  scalar  in  time.  The  stochastic  difference  equation  introduces  an  exchange  rate  7„ 
which  we  model  in  a first  step  as  a deterministic  function.  In  a second  step,  we  generalize 
7„  as  a stochastic  variable  taking  fluctuations  in  the  inhomogeneous  environment  into 
account.  In  Section  3 we  solve  the  non-linear  integral  equation  numerically  and  analyze 
the  influence  of  the  different  parameters  on  the  decay  rate.  The  paper  finishes  with  a 
conclusion. 
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Figure  1.  Mixing  problem  along  Lagrangian  path  way. 


2.  Mixing  model 


2.1.  Problem  statement 


The  evolution  of  a passive  scalar  c = c(t,  x)  is  governed  by  the  stochastic  PDE 


- + v(t,x)yc  = DV^c,  (2.1) 

where  v{t,  x)  is  the  random  velocity  field  and  D is  the  molecular  (or  thermal)  diffusivity. 
The  classical  problem  is  to  derive  a closed  equation  for  the  Euler  one-point  probability 
density  function  (pdf)  p = p{c\t,x)  = {6{c  — c{t,x))).  A detailed  discussion  of  this, 
still  unsolved,  problem  can  be  found  in  Pope  (2000).  In  the  present  paper  we  consider 
the  mixing  problem  in  the  Lagrangian  framework  by  introducing  a passive  scalar  c{t)  = 
c{t,x{t))  of  a particle  volume  moving  with  the  velocity  v{t,x{t))  (see  Fig.  1). 


2.2.  Stochastic  difference  equation 

The  equation  for  the  scalar  c{t)  = c{t,x{t))  of  the  particle  volume  moving  with  the 
random  velocity  v(t,  x{t))  can  be  approximated  by  the  stochastic  equation 

^ = -7(0(c-c(0).  (2.2) 

where  the  exchange  rate  j{t)  and  the  ambient  concentration  c{t)  are  random  processes. 
For  the  constant  values  of  7 and  c,  we  obtain  the  lEM  model  (1.1)  with  the  exponential 
decay: 

= exp{-7i}  . (2.3) 

The  key  feature  of  the  present  model  is  that  c{t)  is  a random  process  and  functional  of 
c{t)  itself.  The  crucial  assumption  that  relates  Eq.  (2.2)  to  the  law  of  large  numbers  is 
that  the  mean  value  {'j{t))  behaves  as  t~^  for  large  t.  One  can  show  that  for  a Gaussian 
scalar  field  c = c{t,x),  the  mean  value  (7(f))  is  proportional  to  the  ratio  where  x 

is  the  mean  scalar  dissipation  rate,  \ = 2jO((Vc)^)  (Pope  2000).  It  follows  from 


M = 

dt 


(2.4) 
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that  if  at"  then  x decays  as  t and  therefore  (7(f))  at  The  appearance 

of  the  power  law  decay  can  be  understood  if  we  assume  that  7 = A(to  + 1)“^ 


dc  _ A 
dt  to  T t 

The  solution  to  this  simplified  equation  is 


c). 


c(t)  - c 
c(0)  — c 


(2.5) 


(2.6) 


Since  our  aim  here  is  to  relate  the  mixing  problem  to  the  law  of  large  numbers,  it  is  more 
convenient  to  rewrite  Eq.  (2.2)  as  a stochastic  difference  equation 


Cn+I  - C„  = -7n(Cn -Cn+i)  , 0 < c„  < 1 , fom  = 1, 2, 3, . . . (2.7) 


where  7„  and  c„,  are  assumed  to  be  sequences  of  mutually  independent  random  variables 


with  the  densities 

w(n,7)  = ^-P{7n  <7}  , 

and  the  first  moments 

1 

(7n)  = / 7^(n,7)d7  = — — 
J 1 4-  n 

0 


4>{n,c)  = 

= 4p{5.<s} 

(2.8) 

(dn)  = 

1 

f c cj>{n,'c)dc  = fj, . 

(2.9) 

0 


Here  the  mean  exchange  parameter  (7„)  has  been  chosen  in  such  a way  that  it  is  equal 
to  A/2  at  time  n = 1.  The  parameter  A can  be  regarded  as  a measure  of  the  mixing 
intensity.  Since  0 < c„  < 1,  7„  obeys  the  inequality  0 < 7n  < 1.  By  introducing  the 
deviations  from  the  mean  //  for  the  concentration  of  the  particle  and  its  surrounding, 
respectively. 


— Cn  M ) (2.10) 

Cn  =Cn~IJ,  (2.11) 

we  can  rewrite  Eq.  (2.7)  as 

^n+l  “ ^ni'^Ti  ^n+l)  j ^ ^1  M , 77  = 1,  2,  3,  . . . , (2.12) 

where  is  a sequence  of  zero  mean,  independent  random  variables  with  the  density 
4>{n,  and  — /z  < $ < 1 — /r. 


2.3.  Law  of  large  numbers  and  forward  Kolmogorov  equation 
To  illustrate  the  connection  between  the  mixing  problem  and  the  law  of  large  numbers, 
consider  the  case  when  the  sequence  is  stationary,  that  is  ip{^)  is  independent  of  n, 
and  the  exchange  parameter  7„  is  a deterministic  sequence  of  the  form 


7n  = 


1 

1 + n 


(2.13) 


with  A = 1.  If  we  assume  uj  - ^i,  it  follows  from  Eq.  (2.12)  that  = (^i  + ^2)/2,ti3  = 
(^1  + ^2  + ^s)/3  and  so  on.  Therefore,  the  solution  of  the  equation  (2.12)  can  be  written 
as  a sample  mean 


1 " 
‘^n  — ~ ^ 

n 


fc=l 


(2.14) 
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which  tends  to  the  mean  /x  as  n — > oo  while  the  variance  (u^)  oc  n~^. 

The  advantage  of  having  an  equation  (2.12)  is  that  we  can  easily  derive  an  equation 
for  the  probability  density  function 


p{n,u)  = < w}  , 


n = 1, 2, 3, . . . . 


(2.15) 


It  follows  from  Eqs.  (2.12)  and  (2.13)  that  u„  is  a discrete  Markov  processes  and  the  pdf 
for  this  process  satisfies  the  forward  Kolmogorov  equation  (Feller  1966) 


p{n  + l,u) 


n + 1 y f 


{n  + l)u-^ 


n = l,2,3,....  (2.16) 


Here  we  used  the  inverse  equation 

Ufi  ~ n [(1  . 


(2.17) 


The  main  properties  of  the  solution  of  the  Kolmogorov  equation  (2.16)  with  the  arbitrary 
initial  condition  p(l,tt)  are 


p(n,u)  ^ S{u)  as 


oo 


and 


1-/J 

I 


u^p{n,  u)du  — » n 


-1 


as 


oo 


(2.18) 


(2.19) 


-M 


The  asymptotic  behavior  of  p{n,  u)  for  large  n is  quite  universal.  It  can  be  written  as 

p{n,u)  oc  -\/nexp  {— n5(u)}  for  n 3>  1 , (2.20) 

where  the  function  S{u)  depends  on  the  particular  choice  of  the  density  (p(^)  for  the 
random  sequence  ^k-  One  can  show  that  S{u)  obeys  the  equation  (Knessl  et  al.  1985) 

i-n 

u^-5  = lnM^^j,  M{x)=  J exp{x^}<p(0d^ . (2.21) 

-M 

2.4.  General  case:  random  exchange  rate  and  non-linear  equation  for  the  pdf 

To  account  for  the  entire  spectrum  of  time  scales  we  now  assume  that  the  exchange 
parameter  7„  is  a stochastic  variable.  It  is  convenient  to  write 

A 


(2.22) 


where  £„  is  the  sequence  of  zero  mean,  independent  random  variables.  Since  0 < 7„  < 1, 
it  follows  from  Eq.  (2.22)  that 


The  decay  of  the  variance  of  the  exchange  rate,  <r^{n),  can  be  determined  as 


(2.23) 


= ((7n  - (7n))^)  = 


(1  +n)2  ■ 


(2.24) 
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When  both  and  are  random  variables,  we  have  to  specify  the  joint  probability 
density  for  e and  If  we  denote  it  by  ^(n,  e,^),  then  the  forward  Kolmogorov  equation 
for  p{n,  u)  takes  the  form 


p(n  + l,.)=/  / + 


-n  -1 


+ l)u  — A(1  + e)^ 
+ n — A(1  + e) 


^{n,£,^)d£dC 


(2.25) 

This  equation  follows  from  the  stochastic  difference  equation  for  the  Markov  jump  process 

u„+i  =Un-  r-^(l  + £n)(w„  - ^n+i)  ~ M < '“n  < 1 “ (2-26) 

1 + n 

and  its  inverse  function 


Un  — 


(n  + l)ttra+l  ~ A(1  + g)^n+l 
1 -|-  n — A(1  d“  s) 


(2.27) 


Since  7„  is  determined  by  the  random  scalar  dissipation  rate,  it  is  natural  to  assume 
that  7„  and  are  independent.  Moreover,  the  density  for  must  be  related  to  p{n,u). 
The  simplest  choice  would  be  to  assume  p{n,^)  which  gives  ^'(n,  £,^)  = p{n,^)ip{n,£) 
and  the  generalized  non-linear  integral  equation  for  p(n,  u)  can  be  written  as 

p(n + !,„)=/  I — 

-fi  -1 

(2.28) 

for  n = 1, 2, 3, It  follows  from  Eq.  (2.24)  that  the  variance  of  the  exchange  rate  <T^(n) 

decays  as  n~^.  In  the  subsequent  sections  we  assume  for  simplicity  that  the  variance 
erf  = (ef ) is  constant. 

We  expect  that  the  variance  of  a passive  scalar  erf  (n)  = ((c„  — /i)^)  = (uf ) behaves  as 
n~“  for  large  n.  Taking  the  analytical  solution  of  the  simplified  model  (2.5)  into  account, 
we  anticipate  a strong  dependence  of  a on  the  mixing  intensity  A. 


^n. 


(n+l)u  - A(1  -|-£)g 
1 + n — A(1  + £) 


'^p{n,^)ip{n,£)d£d^ 


3.  Results 

The  integral  equation  (2.28)  describes  the  relaxation  of  an  arbitrary  initial  density 
distribution  of  a conserved  scalar  c to  a ^-distribution  for  large  n.  Because  of  the  non- 
linearity in  c and  c,  the  equation  needs  to  be  solved  numerically  for  non-trivial  initial 
conditions. 

The  trivial  case  A = 0 describes  the  motion  of  an  inert  particle  {e.g.  dye)  in  its  surround 
and  Eq.  (2.12)  has  the  simple  solution  Un  = ui  for  all  n.  However,  if  A is  close  to  its  upper 
limit,  the  Lagrangian  particle  volume  experiences  strong  interactions  with  its  surrounding 
and  intuitively,  it  can  be  expected  that  the  mixing  intensity  A has  a distinct  influence  on 
the  decay  of  erf  (see  solution  (2.6)  to  simplified  equation). 

The  stochastic  variable  £ influences  the  exchange  rate  7„  and  the  effect  of  its  dis- 
tribution on  the  decay  exponent  a is  desirable  to  understand.  Formally,  we  can  ex- 
press a as  a non-explicit  available  function  of  the  form  a = a(X,tp{n,£))  where  A = 
A(^u/^c)Pe  or  Re).  In  the  subsequent  sections  the  behavior  of  Eq.  (2.28)  for  different 
parameters  and  initial  conditions  is  analyzed. 
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3.1.  Numerical  issues 

In  order  to  advance  the  density  p(n  + l,w)  in  time,  Eq.  (2.28)  must  be  integrated  over 
the  and  £-space  which  results  in  the  overall-complexity  0(iV„  x x N^)  for  each 
time  step.  Here,  denotes  the  number  of  grid  points  in  the  ^-direction.  This  makes  the 
computation  rather  time-consuming.  However,  a computer  code  can  be  easily  parallelized 
along  the  u-coordinate. 

The  initial  density  of  u relaxes  to  a (S-distribution  with  steep  gradients  around  u = 0. 
An  adaptive  grid  using  the  equidistribution  principle  is  employed  in  u-  and  ^-space  in 
order  to  resolve  the  shape  of  the  pdf  properly.  In  the  equidistribution  principle  (Liseikin 
1999),  an  initially  given  number  of  grid  points  is  distributed  in  such  a way  that  the 
condition  w(u)du  — (Nu  - 1)“^  w{u)du  is  satisfied.  The  weight  function  w{u) 
is  a positive  function  measuring  the  variation  of  the  solution. 

3.2.  Parameters  and  initial  conditions 

Prom  experiments  and  DNS  studies  it  is  shown  that  the  density  distribution  of  a passive 
scalar  can  be  reasonably  approximated  by  a ^^-distribution.  The  ;8-distribution  has  non- 
zero probabihty  in  the  interval  [un,in,Umax]  and  has  the  form 

~ «min)^~““'’(w  - Umin)“"H^max  ~ . (3.1) 

It  is  fully  described  by  two  parameters,  which  are  functions  of  (u)  and 

1 , (3.2) 

(3.3) 

The  /3-distribution  adopts  a wide  range  of  shapes:  for  the  maximum  variance  cr„  = //(I— /i) 
and  zero  mean  a double  J-distribution  at  Wmin  and  u„iax  is  approached  and  for  small 
a Gaussian  distribution  around  (u)  is  obtained.  Throughout  the  following  simulations, 
this  distribution  is  used  as  initial  condition  for  p(l,u). 

In  order  to  solve  Eq.  (2.28)  we  have  to  provide  an  appropriate  density  distribution 
for  the  random  fluctuation  e.  The  time-dependent  distribution  of  e is  certainly  problem- 
dependent  and  unknown  in  the  present  consideration.  Therefore  we  assume  for  simplicity 
that  £ has  a statistically  stationary  distribution.  A reasonable  choice  would  be  to  assume 
that  e is  distributed  following  Eq.  (3.1).  A particulair  case  is  obtained  if  V’(e)  has  a 
(5-distribution.  For  the  following  test  cases  we  use  both  distributions  as  prescribed  prob- 
abilities. The  constant  parameters  for  all  simulations  of  Eq.  (2.28)  are  summarized  in 
Table  1. 

3.3.  Influence  of  initial  distribution  inp{l,u) 

In  the  first  numerical  experiment  the  dependence  of  the  decay  exponent  a on  the  initial 
conditions  in  u and  e is  investigated.  Therefore,  prescribed  distributions  with  different 
initial  variance  for  a2(l)  and  are  used  (see  Table  2 and  Fig.  2).  For  this  simulation 
A is  chosen  to  be  unity  and  e is  in  the  range  — 1 < £ < (1  -|-  n)/A  — 1.  The  maximum 
possible  value  for  £ is  then  found  for  n = 1. 

The  decay  rate  of  and  the  decay  exponent  a are  evaluated  and  plotted  in  Fig.  3. 
In  the  left  column,  the  four  different  combinations  with  initial  variance  cr^{l)  = 0.20  are 
plotted,  the  middle  column  shows  the  decay  rate  for  the  cases  with  cr„(l)  = 0.15.  In  the 


(u)  r ((^)  ^niin)(^max  (^}) 


{u)  - Umin  [((^)  -Mmin)(^max  ~ (^)) 

■^max 
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Parameter 

Description 

Value 

nEnd 

number  of  time  steps 

100,000 

P- 

mean  value 

0.5 

iVu 

dimension  of  u 

101 

Ni 

dimension  of  $ 

101 

Nc 

dimension  of  £ 

101 

Table  1.  Fixed  parameters  used  for  the  numerical  simulation. 


Run  p(l,u)  =iIj{£) 


BBl 

/3(0.20) 

/3(0.500) 

BB2 

/3(0.20) 

/3(0.250) 

BB3 

/3(0.20) 

/3(0.125) 

BB4 

/3(0.15) 

/3(0.500) 

BBS 

/3(0.15) 

/3(0.250) 

BB6 

^(0.15) 

^(0.125) 

Run  p(l,u)  ii>{\ , £)  = i){£) 


BBT 

/3(0.05) 

yS(0.500) 

BBS 

;S(0.05) 

/3(0.250) 

BB9 

/3(0.05) 

/?(0.125) 

BDl 

/3(0.20) 

5(e) 

BD2 

/3(0.15) 

5(e) 

BD3 

/3(0.05) 

5(e) 

Table  2.  Combinations  of  initial  distributions  used  for  the  parameter  study:  /3(cT{)  - 
^-distribution  with  prescribed  variance  er^;  5(C)  - Dirac  5-function  for  C- 


Figure  2.  Initial  conditions  for  p(l,u)  and  ^(l,e)  = — 1-0- 


right  column,  the  evolution  of  p(n,  u)  is  shown  for  the  mono-modal  initial  /3-distribution 
with  o"u(l)  = 0.05.  Prom  the  bottom  row  of  Fig.  3 it  can  be  seen  that  p(n,  u)  experiences 
strong  dynamics  over  the  initial  interval,  say  for  n < 10.  Over  this  interval  the  variance 
does  not  follow  the  power  law  decay.  Furthermore  it  can  be  deduced  that  the  stochastic 
perturbation  e reduces  the  decay  rate  of  if  is  close  to  cr^ ^nax  curves). 

The  influence  of  the  initial  distribution  of  e on  the  evolution  of  p[n,  u)  and  consequently 
on  (T^{n)  vanish  for  large  n.  This  is  most  pronounced  for  the  cases  BB7-BB9  and  BD3 
(right  column).  Here  the  decay  exponents  for  all  of  these  cases  collapse  to  a single  curve, 
independent  of  ■0(e). 


Rim  p(l.w)  V’(s) 

Dynl  [<5(/x  + w)  + J(1  — p — ii)]/2  (5(e) 

Dyn2  y9(0.15)  ^(0.50) 


Table  3.  Combinations  used  to  study  dynamic  behavior  of  p{n,  u)  for  A = 0.75. 


Prom  this  simulation  it  can  be  concluded  that  the  initially  assumed  distribution  for  e 
has  only  marginal  influence  on  the  decay  exponent  a for  large  n.  During  the  simulation 
the  value  for  a never  reaches  a steady  state.  The  final  convergence  rate  for  n = 100, 000 
is  |do:/dn|  < 1.0  x 10“®  so  that  a is  in  the  range  1.50  < a < 1.65  for  all  test  cases  used 
in  the  present  section. 

3.4.  Dynamic  behavior  ofp{n,u) 

In  the  present  section  we  are  interested  in  the  influence  of  ip{e)  on  the  dynamics  of 
p[n,v).  The  transitional  behavior  of  can  be  retarded  if  the  mixing  intensity  A is 
small.  This  effect  was  already  discussed  previously.  In  the  present  section  the  mixing 
intensity  A = 0.75  is  used. 

In  order  to  obtain  insight  into  the  transitional  behavior  we  use  a double  (5-distribution 
and  a bimodal  /3-distribution  as  initial  condition  for  u,  p(l,u)  = [5{jjL+u)+5{\— p—u)]/2 
and  p{l,u)  = /3(0.15).  For  ip{e)  = 6{e)  a simple  analytical  solution  to  Eq.  (2.28)  can  be 
found  for  the  first  time  step  n = 1.  For  clarity  the  different  combinations  used  here  are 
summarized  in  Table  3.  At  n = 1 the  initial  distribution  for  the  case  Djml  splits  up  into 


Figure  4.  Transitional  behavior  of  p(n,u)  for  different  initial  conditions  and  mixing  intensity. 


Run 

A 

p(l,u) 

V>(e)  [t 

-min  j £max] 

LAI 

1.50 

/3(0.15) 

)3(0.20) 

[-1,1/3] 

LA2 

1.00 

/3(0.15) 

)9(0.20) 

[-1,1] 

LAS 

0.75 

/3(0.15) 

;3(0.20) 

[-1,5/3] 

LA4 

0.50 

/3(0.15) 

;9(0.20) 

-1,3] 

LAS 

0.25 

;d(0.15) 

;3(0.20) 

-1,7] 

Table  4.  Combinations  used  to  study  the  influence  of  the  mixing  intensity  A on  the  evolution 

of  the  variance. 


four  peaks  at  positions  u = A(^  — 7)/2  +7  with  7 = {—p,  1 — m}  and  f = {— Mi  1 ~ m}-  In 
the  succeeding  time  those  peaks  move  toward  the  center  u = 0 whereas  the  tails  flatten. 

An  entirely  different  dynamic  behavior  is  obtained  for  the  case  Dyn2.  The  initially 
smooth  density  function  p(l,u)  transits  rapidly  to  a mono-modal  distribution  and  con- 
verges to  the  5-pdf.  With  the  results  of  this  and  the  previous  section  we  can  summarize 
that  the  prescribed  distribution  of  e influences  the  transitional  dynamics  of  p{n,  u),  how- 
ever its  effect  on  the  decay  exponent  a for  large  n is  insignificant. 

3.5.  Mixing  intensity  A 

From  experiments  and  DNS  studies  it  is  well-known  that  the  decay  rate  strongly  depends 
on  the  initial  ratio  of  the  velocity  and  scalar  length-scales,  tu  and  tc-  The  free  param- 
eter A is  introduced  in  the  exchange  rate  (2.22)  in  order  to  link  the  present  model  to 
physical  processes.  It  needs  to  be  emphasized  that  we  did  not  attempt  to  reproduce  any 
experimental  or  numerical  data,  amply  and  in  great  diversity  available  in  the  literature. 
Presently,  we  Eire  only  interested  in  understanding  the  general  behavior  of  Eq.  (2.28).  It 
was  concluded  previously  that  the  presumed  form  of  tj}{£)  has  only  insignificant  influence 
on  the  decay  rate.  However,  the  dynamic  behavior  of  p{n,  u)  in  the  transitional  range 
is  influenced  by  the  stochastic  variable  e.  In  anticipation  of  the  overall  weak  influence 
of  ^(s)  on  the  decay  exponent  we  chose  arbitrarily  ^(e)  = /3(0.20)  whereas  ^(e)  has 
non-zero  probability  only  in  the  interval  — l<e<2/A  — 1.  The  parameters  used  for  the 
five  different  test  cases  are  summarized  in  Table  4. 

Figure  5 shows  the  distinct  dependence  of  the  mixing  intensity  on  the  decay  of  aj. 


Stochastic  mixing  model  with  power  law  decay  of  variance 
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Figure  5.  Decay  of  cr^  for  different  mixing  intensities  A. 


Experiments  on  turbulent  mixing  show  that  the  variance  decreases  faster  with  in- 
creasing ratio  £u/^c>  until  about  ^u/4  ~ 5.  For  ratios  greater  than  that  value  the  depen- 
dence is  negligible.  This  tendency  can  be  resembled  by  changing  the  mixing  intensity  A. 
The  power  law  exponent  a is  approximately  1.77  for  high  mixing  intensity  and  decreases 
to  about  0.62  for  A = 0.25  which  corresponds  to  an  indolent  mixing  process.  From  the 
present  findings  it  is  difficult  to  specify  the  exact  functional  relation  between  A and  initial 
conditions  used  in  experiments.  Therefore  we  express  the  mixing  intensity  formally  as 
A = X{£u/^c,  or  Pe). 


4.  Conclusion 

In  the  present  paper  we  presented  a robust  mixing  model  based  on  the  law  of  large 
numbers.  The  exchange  rate  decays  as  t~^  and  links  the  mixing  model  to  the  law  of  large 
numbers.  By  adding  stochastic  fluctuations  to  the  exchange  rate,  we  extend  the  model 
to  the  general  case,  taking  the  fluctuations  of  the  surround  into  account. 

The  sensitivity  of  the  model  is  analyzed  for  different  initial  conditions  and  parameters. 
It  was  shown  that  the  initial  distribution  for  e and  u has  only  insignificant  influence  on 
the  decay  rate  of  the  variance.  This  is  to  be  expected  because  the  initial  distribution  of 
the  scalar  is  almost  never  measured  in  experiments,  assuming  implicitly  that  either  the 
scalar  is  distributed  initially  as  a double  5-function  or  generally  disregard  the  dependence 
of  the  decay  of  on  the  initial  density  distribution. 

It  was  shown  that  with  increasing  mixing  intensity  the  variance  of  the  passive  scalar 
decreases  faster.  The  connection  of  our  model  to  physical  mixing  processes  is  provided 
through  A,  which  can  account  for  initial  length-scale  dependences  £u/'^c  and  other  physical 
effects. 
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